arXiv:1501.02411vl [cs.LG] 11 Jan 2015 


1 


A Gaussian Particle Filter Approach for 
Sensors to Track Multiple Moving Targets 
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Abstract —In a variety of problems, the number and state 
of multiple moving targets are unknown and are subject 
to be inferred from their measurements obtained by a 
sensor with limited sensing ability. This type of problems 
is raised in a variety of applications, including monitoring 
of endangered species, cleaning, and surveillance. Particle 
filters are widely used to estimate target state from its prior 
information and its measurements that recently become 
available, especially for the cases when the measurement 
model and the prior distribution of state of interest are 
non-Gaussian. However, the problem of estimating number 
of total targets and their state becomes intractable when 
the number of total targets and the measurement-target 
association are unknown. This paper presents a novel 
Gaussian particle filter technique that combines Kalman 
filter and particle filter for estimating the number and state 
of total targets based on the measurement obtained online. 
The estimation is represented by a set of weighted particles, 
different from classical particle filter, where each particle 
is a Gaussian distribution instead of a point mass. 

Index Terms —Kalman Filter, Particle Filter, Gaussian 
Particle Filter, Multiple Target Tracking 

1. INTRODUCTION 

The problem of tracking and monitoring targets us¬ 
ing position-fixed sensors is relevant to a variety of 
applications, including monitoring of moving targets 
using cameras [], tracking anomalies in manufacturing 
plants Q, and tracking of endangered species 0-0. 
The position-fixed sensor is deployed to measure targets 
based on limited information that only becomes available 
when the target enters the sensor’s field-of-view (FOV) 
or visibility region 0 - The sensor’s FOV is defined as 
a compact subset of the region of interest, in which 
the sensor can obtain measurements from the targets. 
In many such sensor applications, filter techniques are 
often required to estimate unknown variables of interest, 
for examples, target number and target state. 

When the noise in the measurement model is an 
additive Gaussian distribution, the target state can be 
estimated from frequent observation sequence using a 
Kalman filter 0 - This approach is well suited to long- 
range high-accuracy sensors, such as radars, and to 
moving targets with a known dynamical model and 
initial conditions. However, most of these underlying 
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assumptions are violated in modern applications, because 
the targets’ motion models are unknown, and, possibly, 
random and nonlinear. Also, due to the use of low-cost 
passive sensors, measurement errors and noise may be 
non-additive and non-Gaussian. An extended Kalman 
filter (EKF) can be used when the system dynamics 
are nonlinear, but can be linearized about nominal 
operating conditions Q. An unscented Kalman filter 
(UKF) method, based on the unscented transformation 
(UT) method, can be applied to compute the mean and 
covariance of a function up to the second order of the 
Taylor expansion jSj. 

However, the efficiency of these filters decreases sig¬ 
nificantly when the system dynamics are highly nonlin¬ 
ear or unknown, and when the measurement noise are 
non-Gaussian. Recently, a non-parametric method based 
on condensation and Monte Carlo simulation, known as 
a particle filter, has been proposed for tracking multiple 
targets exhibiting nonlinear dynamics and non-Gaussian 
random effects ©• Particle filters are well suited to 
modem surveillance systems because they can be applied 
to Bayesian models in which the hidden variables are 
connected by a Markov chain in discrete time. In the 
classical particle filter method, a weighted set of particles 
or point masses are used to represent the probability 
density function (PDF) of the target state by means of 
a superposition of weighted Dirac delta functions. At 
each iteration of the particle filter, particles representing 
possible target state values are sampled from a proposal 
distribution |T0| . The weight associated with each par¬ 
ticle is then obtained from the target-state likelihood 
function, and from the prior estimation of the target 
state PDF. When the effective particle size is smaller 
than a predefined threshold, a re-sampling technique 
can be implemented (TT). One disadvantage of classical 
particle-filtering techniques is that the target-state transi¬ 
tion function is used as the importance density function 
to sample particles, without taking new observations into 
account (ID Recently, an particle filter with Mixture 
Gaussian representation was proposed by author for 
monitoring maneuvering targets GD’ where the particles 
are sampled based on the supporting intervals of the 
target-state likelihood function and the prior estimation 
function of the target state. In this case, the supporting 
interval of a distribution is defined as the 90% confidence 
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interval 0. The weight for each particle is obtained 
by considering the likelihood function and the transition 
function simultaneously. Then, the weighted expectation 
maximization (EM) algorithm is implemented to use the 
sampled weighted particles to generate a normal mixture 
model of the distribution. 

Kreucher proposed joint multitarget probability den¬ 
sity (JMPD) p5| to estimate the number of total targets 
in workspace and their state, where targets are moving. 
By using JMPD, the data association problem is avoided, 
however, the JMPD results in a joint state space, the 
dimension of which is the dimension of a target state 
times number of total targets. Since the number of 
total targets is unknown, the joint space size remain 
unavailable. To overcome this problem, it is assumed the 
number of total targets has a maximum value. Therefore, 
when the maximum number of targets is large, the joint 
state space becomes intractable. 

Inspired by | [T6| , this paper presents a novel filter tech¬ 
nique which combines Kalman filter and particle filter 
for estimating the number and state of total targets based 
on the measurement obtained online. The estimation is 
represented by a set of weighted particles, different from 
classical particle filter, where each particle is a Gaussian 
instead of a point mass. The weight of each particle 
represents the probability of existing a target, while its 
gaussian indicates the state distribution for this target. 
More importantly, the update of particles is different 
from classical particle filter. For each particle, the gaus¬ 
sian parameters are updated based using Kalman filter 
given a measurement. To overcome the data association 
problem, in this paper, when one particle is updated, 
the other particles are considered as the measurement 
condition, which will be explained in Section IV The 
novel Gaussian particle filter technique requires less 
particles than classical particle filters, and can solve 
multiple target estimation problem without increasing the 
state space dimensions. 

The paper is organized as follows. Section |IJ de¬ 
scribes the multiple targets estimation problem formu¬ 
lation and assumptions. The background on the particle 
filter and Kalman filter is reviewed in Section [nil Section 


IV presents the Gaussian Particle filter technique. The 


method is demonstrated through numerical simulations 
and results, presented in Section fV] C onclusions and 
future work are described in Section IvTI 


II. Problem Formulation 

There N targets are moving in a two dimensional 
workspace denoted as W where N denotes the unknown 
number of total targets. For simplicity, there is zero 
obstacle in the workspace. The goal of the sensor is to 
obtain the state estimation for all the targets, denoted 
as X/c, and target number estimation, denoted as Tk, at 


time step k. The states for total targets at k is denoted 
as X/c = ,x^] that has N state vectors. 

The estimation of total target state at k is denoted as 
X^ = [x^,x^, • • • ,x^ ]. Let denote the estimation 
of total target number at The ith target is modeled 
as 

Xl = F/eX*^_i +IZ/,, (1) 

where 

iZfc-7V(0,Qfc). (2) 

Furthermore, Fj^ and Qk are assumed known. 

In standard estimation theory, a sensor that obtains a 
vector of measurements G in order to estimate an 
unknown state vector set X^ G at time k is modeled 
as, 

= h(X^A^), (3) 

where h : ^ is a deterministic vector 

function that is possibly nonlinear, the random vector 
\k ^ represents the sensor characteristics, such as 
sensor action, mode GI), environmental conditions, and 
sensor noise or measurement errors. In this paper, the 
sensor is modeled as 

1 ^ 

Zk =+ ujk- ( 4 ) 

^ i=l 

It is further assumed that the whole workspace is visible 
to a position fixed sensor (not shown). 

III. Background 
A. Particle Filter Methods 

The particle filter is a recursive model estimation 
method based on sequential Monte Carlo Simulations. 
Because of their recursive nature, particle filters are 
easily applicable to online data processing and variable 
inference. More importantly, it is applicable to nonlinear 
system dynamics with non-Gaussian noises. The PDF 
functions are represented with properly weighted and 
relocated point-mass, known as particles. These particles 
are sampled from an importance density that is crucial 
to the particle filter algorithm and is also referred to 
as a proposal distribution. Let denote 

the weighted particles that are used to approximate the 
posterior PDF /(x^ | Z'^) for the jth target at where 
Zj = {z^,..., z^} denotes the set of all measurements 
obtained by sensor i, from target j, up to Then, the 
posterior probability density function of the target state, 
given the measurement at can be modeled as, 

N N 

p=l p=l 

where ^ is non-negative and 6 is the Dirac delta 
function. The techniques always consist of the recursive 
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propagation of the particles and the particle weights. In 
each iteration, the particles are sampled from the 
importance density g(x). Then, weight is updated 
for each particle by 




( 6 ) 


where p{'Kj p) oc /(x^^ | Z^). Additionally, the weights 
are normalized at the end of each iteration. 

One common drawback of particle filters is the de¬ 
generacy phenomenon 0, i.e., the variance of particle 
weights accumulates along iterations. A common way 
to evaluate the degeneracy phenomenon is the effective 
sample size GD, obtained by, 

N,. = -- (7) 


Ep=i(^ 


3,P 


where p = 1, 2,..., are the normalized weights. 
In general, a re-sampling procedure is taken when Ne < 
Ns, where Ng is a predefined threshold, and is usually 
set as Y- Let denote the particle set 

that needs to be re-sampled, and let 
denote the particle set after re-sampling. The main idea 
of this re-sampling procedure is to eliminate the parti¬ 
cles having low weights by re-sampling 
from with the probability of p(x^* = 

x^ ^) = At the end of the resampling procedure, 
Wj^p^p = 1, 2,..., are set as 1/N. 


B. Kalman Filter Methods 

The well known Kalman filter is also a recursive 
method to estimate system/target state based on a mea¬ 
surement sequence, minimizing the estimation uncer¬ 
tainty. The measurement of the system state with an addi¬ 
tive Gaussian noise are given by the sensor. Then, in each 
iteration, the Kalman filter consists of two processes: i) 
it predicts the system state and their uncertainties; ii) 
it updates the system state and uncertainties with the 
measurement that newly becomes available. The system 
dynamics is given as 

X/c = F/cX/e_i + BkUk + V'k (8) 

where subscript k and k — 1 denote the current and 

previous time index, while is the system discrete 
transition matrix, and and Uk are the control matrix 
and control input, i^k is the white noise, defined as 

Uk-N{0Mk) (9) 

where S/c is the covariance. At kth time step, an mea¬ 
surement of the system true state x/^ is made by a sensor, 
is given by 

Zk = H/eX/e + U3k (10) 


where FLk is a mapping from system state space to 
measurement space, and white noise Q/c is defined as 

u;k^N{0,Rk) (11) 

It is assumed that the noise uJk and u at each time step 
are independent. 

Let x/c denote the predicted state estimation given 
x/c_i, where x/^-i is the updated estimation of system 
state at k — 1 time step. Furthermore, let J^k denote the 
predicted covariance given where ^k-i is the 

updated estimation covariance. Then, in the predicting 
step, 

x/e = F/eX/c-i + BkUk (12) 

^k = F/cE/c-iF^ + Qk (13) 

In the updating step, the measurement Zk is used, to¬ 
gether with above predicted state and covariance, to up¬ 
date the state and covariance. The residual, yk between 


measurement and predicted state is given by 

yk = Zk- H/eX/c (14) 

The innovation covariance Sk is given by 

Sfc = HfcSfcHj + Rfe (15) 

Then, the optimal Kalman gain is calculated as 

Kfe = tkHlSp (16) 

Then, the state and covariance can be updated by 

x/e = x/e + Kkjk (17) 

±k = {l-KkBk)tk (18) 

IV. Methodology 


In this paper, a novel Gaussian particle filter technique 
follows the main idea of particle filter for estimating 
the number and state of total targets based on the 
measurement obtained online. Different from classical 
particle filter, each particle here is a gaussian instead of 
a point mass. The estimation for number of total targets 
and their state is presented by a set of weighted particles. 
The ith particle at time k is denoted as 

= (19) 

where Wi is the probability of existing a target having a 
state distribution as A^(/L4, A/’(xi|/L4, 5]^)). By this particle 
definition, the dimensions of the system state remains the 
same as the dimensions of each individual target. When 
these particles are available, the estimated number of 
total targets can be given as T = Wi, where Np is 
the particle number. 

Notice that the particle representation is different from 
classical particle filter, where each particle represents a 
possible value of system state. The updating of each 
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particle and total weights are also different from classical 
particle filter. Kalman filter is used to update each 
particle, the weight and the distribution. Notice that since 
the measurement at each time step is conditioned on all 
the targets in the FOV, while in the classical Kalman 
filter method one measurement is associated with one 
target, which means data association problem is avoid. 
Therefore, the Kalman filter is modified to updated 
the particles which are coupled by one measurement, 
and some approximations and assumptions are further 
needed. 


Similar to Kalman filters and particle filters, the al¬ 
gorithm proposed in this paper is a recursive method. 
Assume that at time step k, the measurement Zk is 
available, and the estimation of the system at time step 
k — lis represented by a particle set, denoted SisVk-i = 
P|_^,..., P^^}, where Np is the number of 
all particles. By using the target dynamic function 
Vk-i can be updated to without using the Due to 
limit of FOV, only a few particles may have contribution 
to the measurement. Let Vs denote set including the 
particles lie in the FOV, while let P = VkiVs denote 
the compensation set. Only the particles in Vs are 
updated. Please Note the size of Vs is small. Without 
generality, assume that Vs = {P^, P|,..., P^}, where 
s the number of particles The update of each particle in 
Vs is calculated separately. Without generally, we focus 
on updating P^, Let a boolean set P = [ei, ei,..., e^], 
where G {0,1}. For any E with ej = 1 such that 
> e, where e is a predefined 
threshold, a particle is calculated and denoted as pj, 
Then, the modified Kalman filter is used to give the 
updated gaussian parameters of all particles with = 1. 
According to sensor model the measurement is given 
by 


Zk = 


Etie. 


(EEie, 


^ ejSj+cjfc (20) 


Compare the above function to ( p^ , we have following 
setting 

1 


Hk = Ix 


Rk = 


Zk = Zk- 
1 


^ (EEieO 
EEie. 






( 21 ) 

( 22 ) 

( 23 ) 


Then, by applying Kalman procedure 


Yfc = Zfe - ( 24 ) 

Sfe = + Rfc ( 25 ) 

Kk = ( 26 ) 

Mfe = Mfc + KfeYfe ( 27 ) 

Sfc = (I-KfcHfc)Sfc ( 28 ) 


Its proof can be found in the appendix. 

Once each particle appearing in combination E has 
been updated, the weight Wc is for the particle combi¬ 
nation E can be obtained by 


Wc 


= n(wj)®*(i - Wi)^ X 


(27r)2||s;; 


X exp{-(zfc - He) (Zfc - Me)} ( 29 ) 


where c e Ie, where is the combination index, and 
/j,c) and is given by 

^lc = nkE^^k ( 30 ) 

= ( 31 ) 


Then, insert particle {A/’*(/L6 /c, H/c)} into a set Gc for 
combination c G I^;, the set Gc has a weight Wc. After all 
Gc the combination of E that > e 

is obtained. Weights are updated by 


Wc 


Wc 


( 32 ) 


. Then in each group Gc, the weight of ith particle is 
updated as 


w 


c 



^W 

i 

k 


c 


(33) 


The particles in all set Gc are updated from the same 
particle in the previous set. If two particles in is close 
enough, then they are combined as one particle, the 
weight of which is set as the summation of both weights. 
The distance between and pij is defined as Mahal- 
distance 

(Mi - - IJ-j) ( 34 ) 


and its covariance is updated as 

S = (35) 

V. Simulation and Results 

As shown in figure N targets, represented by 
blue dots, are moving in the 2 dimensional workspace. 
The whole workspace is visible to a position fixed 
sensor (not shown) and the workspace is discretized into 
12 X 12 cells. Each cell represents a 1 x 1 rectangular 
area. The ith cell, denoted as GiG G C, is defined by 
where C is the cell index set and 
(x'f^yf^) and {xf^, yf^) are up left and down right cor¬ 
ner coordinates of the ith rectangular area respectively. 
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Only M cells can be measured at each time step k, 
and they don’t have to be adjacent. The goal of the 
sensor is to estimate the target states and target number 
at time k. In this paper, information value function 
based a divergence is used to select the best M cells 
to measure at each step The estimation of target 
states and target number at time k is represented by joint 
multitarget probability density(JMPD) and it is updated 
after obtaining new measurements. 



same cell, then the detection probability is Pd{T) = 
and the ith sensor measurement at 
time k can be evaluated by 






Pd{T) zf = 1 
l-Pd(T) 4 = 0 


N 

T = ^(4 > xf) n (4 < xf) 

i=i 

n(4 > 4') n (4 < xi^), c = (39) 

= (40) 


where Xj^yj are two position components of G 
and is the target number. Additionally, operators ”>” 
and ”<” return either 1 if true or 0 if false, while ’’D” is 
the Boolean operator ’’and”. For example, as shown in 
figure 1^, when c = k,T = 2, similarly, when c = j{i), 
T = 1(0). 

A snapshot of simulations is shown in Fig. where 
magenta squares represent positive measurement return 
and blue dots represent the true targets’ positions. The 
simulation results are summarized in Fig. where the 
black curve represents the target state estimation error 
and the red curve represents the target number estimation 
error. As shown in Fig. both errors decreases as more 
measurements become available. 


The target time-discrete state transition function can 
be written as 

= Fxf + wf (36) 


where 


1 r 0 0 
0 10 0 
0 0 1 r 
0 0 0 1 


(37) 


and v^^ is 0 mean Gaussian noise with covariance 
Q =diag(20, 0.2, 20, 0.2 ), a nd r is the time step length, 
andi e {1,2,--- ,T4 (^. 

It is further assumed that i) the sensor can measure 
any cell at time k; ii) the sensor can only measure up to 
M cells at time k. The sensor condition Aj represents 
the signal to noise ratio SNR, currently, it has only one 
possible value, fixed and known. The measurement 
is a discrete variable, then joint PMF can be written as 


/(z^ X^ a 4 =/(z''|X^ T^ A4/(X^ t4/(A4 

(38) 

When measuring a cell, the imager sensor will give 
a Raleigh return, either a 0 (no detection) or a 1 
(detection) governed by detecting probability, denoted as 
Pd, and false alarm probability, denoted as p/. According 
to standard model for threshold detection of Rayleigh 
returns, pf = When T targets are in the 



Fig. 2. Simulation Result 


VI. Conclusion and Future Work 

A Gaussian particle filter that combines Kalman fil¬ 
ter and particle filter is presented in this paper for 
estimating the number and state of total targets based 
on the measurement obtained online. The estimation is 
represented by a set of weighted particles, different from 
classical particle filter, where each particle is a Gaussian 
instead of a point mass. The weight of each particle 
represents the probability of existing a target, while its 
Gaussian indicates the state distribution for this target. 
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Fig. 3. Snapshot of simulations 


This approach is efficient for the problem of estimating 
number of total targets and their state. 


VIL Appendix 

Without losing generality, Vs = ..., 

E = [ei, 02 ,..., Cg], for any particle such that ej = 1, 
its fj^i and 5]^, given Zk and Vs- 

y, = z,- (41) 


By setting = 0, therefore 


n = K( 


Then, 


k\ 




-1 


Z^i=l Z^i=l 

4 = iii + Klvt (45) 




3 

k 


COY{xi-xi) 

COy(xi-(xi^Kiy,)) 

i<-i4 + Ki{ 


cov 


jfEUl Xk(^i 


Ek'k — 


e:=i* 


k^^ 


T.Uei 


COV (I 


e:=i' 


EEie. 


-IKi){xi-xi) 


{I - - 

l^i=l 


Xl 


(42) 


EEi®. 


-mr 


y^IKi Y. 

^i=l ^i=l 


j'^T 


+KiR,Ki. 


(43) 


n 




1 


{Rk + 






p. z_^^k\s^s 

Z^i=i z^i= 


1 


(46) 
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